
rm(list=ls(all=TRUE))
library(ggplot2)
library(gridExtra)
library(car)
library(plm)
library(gam)
library(plyr)
library(doBy)
library(foreign)
library(texreg)
library(xtable)
library(MASS)
library(GGally)
library(zoo)
library(lmtest)
library(plm)
library(lattice)
library(ebal)
library(DataCombine)
library(rgenoud)
library(sandwich)
library(Hmisc)
library(dplyr)
library(tidyr)
library(lubridate)
library(zoo)
require(texreg)
require(lmtest)
require(multiwayvcov)
library(VIM)
library(rpart)

setwd("/users/irismalone/Dropbox/AGDPapers/IS-ISSS Paper/Replication Files - II/Data")

df2=read.csv("replication_extsup_data.csv")
df2$l1noviolsupl = log(df2$l1noviolsup+1)
#####################
#Table 1 Summary of Shared Ties
#####################

table(df2$sharedethnicpower)
table(length(unique(df2$torgid[df2$sharedethnicpower==1])))

table(df2$leftideational)
table(length(unique(df2$torgid[df2$leftideational==1])))

#Religious is either Sunni or Shia
table(df2$sunniideational)
table(length(unique(df2$torgid[df2$sunniideational==1])))

table(df2$shiaideational)
table(length(unique(df2$torgid[df2$shiaideational==1])))


#####################
# Hypothesis 1 Results
#####################

reg.m0 = glm((anystatesup) ~ l1noviolsupl + sharedcoideo, data=df2, family="binomial")
vcov_firm <- cluster.vcov(reg.m0, df2$dyadid)
a = coeftest(reg.m0, vcov_firm)
a

reg.m1 = glm((anystatesup) ~ l1noviolsupl*sharedcoideo, data=df2, family="binomial")
vcov_firm <- cluster.vcov(reg.m1, df2$dyadid)
b = coeftest(reg.m1, vcov_firm)
b

reg.m2 = glm((anystatesup) ~ l1noviolsupl*sharedcoideo  + rivalry + minwbgdppc + mindemest, data=df2, family="binomial")
vcov_firm <- cluster.vcov(reg.m2, df2$dyadid)
c = coeftest(reg.m2, vcov_firm)
c

reg.m3 = glm((anystatesup) ~ l1noviolsupl*sharedcoideo + rivalry + minwbgdppc + mindemest + factor(ccode1), data=df2, family="binomial")
vcov_firm <- cluster.vcov(reg.m3, df2$dyadid)
d = coeftest(reg.m3, vcov_firm)
d


texreg(list(a, b, c, d), omit.coef="ccode1", stars=c(0.1, 0.05, 0.01), caption="Main")
texreg(list(reg.m0, reg.m1, reg.m2, reg.m3), stars=c(0.1, 0.05, 0.01))


htmlreg(list(a, b, c, d), omit.coef="ccode1", stars=c(0.1, 0.05, 0.01), caption="Main Results", file = "table2_mainresults.doc")

htmlreg(list(reg.m0, reg.m1, reg.m2, reg.m3), file = "table2_mainresults_unclustered.doc")

#####################
# predicted probability values
#####################

reg.m1 = glm((anystatesup) ~ l1noviolsupl*sharedcoideo + minwbgdppc + mindemest + rivalry, data=df2, family="binomial")
vcov_firm <- cluster.vcov(reg.m1, df2$dyadid)
b = coeftest(reg.m1, vcov_firm)
b
require(ggplot2)
require(ggeffects)

p_1000 <- ggpredict(reg.m1, c("l1noviolsupl [meansd]", "sharedcoideo"), vcov.fun = "vcovCR", vcov.type = "CR0", 
                    vcov.args = list(cluster = df2$dyadid))
p_1000

p_1000 <- ggpredict(reg.m1, c("l1noviolsupl [0, 0.5, 1, 1.5, 2, 2.5, 3, 3.5]", "sharedcoideo"), vcov.fun = "vcovCR", vcov.type = "CR0", 
                    vcov.args = list(cluster = df2$dyadid))
p_1000

a = plot(p_1000, colors="bw") + 
  theme_bw() + ylab("Predicted Pr(External Support)") + xlab(expression("Ln(Number of Groups)"))  + ggtitle("") + 
  scale_colour_manual(values = c(1, 2), name = "Shared Ties") + 
  scale_linetype_manual(values = c(1, 2), name = "Shared Ties") +  
  theme(plot.title = element_text(hjust = 0.5)) +
  theme(legend.position = c(0.7, 0.8))
a
pdf("prednum.pdf", 4, 3)
plot(p_1000, colors="bw") + 
  theme_bw() + ylab("Predicted Pr(External Support)") + xlab(expression("Ln(Number of Groups)"))  + ggtitle("") + 
  scale_colour_manual(values = c(1, 2), name = "Shared Ties") + 
  scale_linetype_manual(values = c(1, 2), name = "Shared Ties") +  
  theme(plot.title = element_text(hjust = 0.5)) +
  theme(legend.position = c(0.7, 0.8))
dev.off()

